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1 INTRODUCTION 



ABSTRACT 

The axial modes for non-barotropic relativistic rotating neutron stars with uniform 
angular velocity are studied, using the slow-rotation formalism together with the low- 
frequency approximation, first investigated by Kojima. The time independent form 
of the equations leads to a singular eigenvalue problem, which admits a continuous 
spectrum. We show that for I = 2, it is nevertheless also possible to find discrete mode 
solutions (the r-modes). However, under certain conditions related to the equation 
£T) . of state and the compactness of the stellar model, the eigenfrequency lies inside the 

continuous band and the associated velocity perturbation is divergent; hence these 
solutions have to be discarded as being unphysical. We corroborate our results by 
explicitly integrating the time dependent equations. For stellar models admitting a 
physical r-mode solution, it can indeed be excited by arbitrary initial data. For models 
admitting only an unphysical mode solution, the evolutions do not show any tendency 
to oscillate with the respective frequency. For higher values of I it seems that in certain 
cases there are no mode solutions at all. 
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^ , Immediately after the discovery of the r-modes being generically unstable with respect to gravitational-wave emission (An- 
dersson 1998a; Friedman & Morsink 1998), it was suggested that they may cause the newly born neutron stars to spin down 
via the emission of gravitational waves (Lindblom, Owen & Morsink 1998; Andersson, Kokkotas & Schutz 1999). Because of 
their surprisingly fast growth times, r-modes should be able to slow down a hot and rapidly spinning newly born neutron 
star during the first months of its existence. There is also work suggesting that the r-mode instability might be relevant for 
old neutron stars in binary systems. This potential relevance for astrophysics has attracted the interest of both the relativity 
and the astrophysical community on various aspects of this subject. For an exhaustive upto-date review, see for instance 
Andersson & Kokkotas (2001) and Friedman & Lockitch (2001). 

Most of the recent work on r-modes is based on Newtonian calculations under the assumption of slow rotation, and the 
effects of gravitational radiation are incorporated through the quadrupole formula. However, it is clear that for a complete and 
quantitatively correct understanding, one has to use the framework of general relativity. Still the slow-rotation approximation 
is well justified, since the angular velocity of even the fastest spinning known pulsar corresponds to a rotational expansion 
parameter of e = M/B? ~ 0.3. The full set of equations in the slow-rotation limit was first given by Chandrasekhar & 
Ferrari (1991) for the axisymmetric case, and by Kojima (1992) for the general case. 

In the non-rotating case, the perturbation equations are decoupled with respect to the harmonic index / and degenerate 
with respect to the azimuthal index m. Furthermore the oscillation modes can be split into two independent sets, which 
are characterized by their behaviour under parity transformation. The polar (or spheroidal) modes transform according to 
( — 1) ! , whereas the axial (or toroidal) modes according to (— l) i+1 . For a non-rotating perfect fluid star, the only possible fluid 
oscillations are the spheroidal /- andp-modes. For non-barotropic stars, i.e. stars with a temperature gradient or a composition 
gradient, there exists another family of modes, the g-modes, where the main restoring force is gravity. All axial perturbations 
of non-rotating perfect fluid stars have zero frequency, i.e. they represent stationary currents. In the non-barotropic case, this 
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zero- frequency space consists only of the axial r-modes, while for barotropic stars [j, it also includes the polar g-modes, since 
they require a temperature gradient for their existence. 

As the star is set into rotation, the picture changes. In the slow rotation approximation, the m-degeneracy is removed 
and the polar modes with index I are coupled to the axial modes with indices I ± 1 and vice versa. Furthermore, the rotation 
imparts a finite frequency to the zero-frequency perturbations of the non-rotating stars. In non-barotropic stars, those modes, 
whose restoring force is the Coriolis force, all have axial parity. However, as has been first pointed out in the Newtonian 
framework by Lockitch & Friedman (1999), for barotropic stars the rotationally restored (inertial) modes are generically 
hybrids, whose limit in the non-rotating case are mixtures of axial and polar perturbations. 

If one focuses only on the r-mode, whose frequency is proportional to the star's angular velocity f2, one can order the 
perturbation variables in powers of il. Kojima (1997, 1998) used this low-frequency approximation, sometimes also called 
slow-motion approximation (Schumaker & Thorne, 1983), to show that the purely axial modes of (non-barotropic) stellar 
models can be described by a single second-order ODE. This eigenvalue problem, however, proves to be singular, since it is 
possible for the highest derivate to vanish at some value of the radial coordinate inside the star. Kojima (1997,1998) then 
argued that this singular structure should give rise to a continuous spectrum. This has been put on a rigorous mathematical 
footing by Beyer & Kokkotas (1999). The appearance of a continuous spectrum can be explained as follows. It is well known 
that in the Newtonian limit, the eigenfrequency of the r-mode for an inertial observer is given by 



A first relativistic correction can be obtained by using the relativistic Cowling approximation, which consists in neglecting 
all metric perturbations. In this case, the only correction comes from the frame dragging ui, which is a function of the radial 
coordinate r, thus leading to an r-dependent oscillation frequency of each fluid layer: 



Instead of a single frequency, there is now a continuous band of frequencies, whose boundaries are determined by the values 
of the frame dragging u)(r) at the centre and the surface of the star. 

However, it has been pointed out (Beyer & Kokkotas 1999) that the existence of the continuous spectrum might be 
just an artefact of the too restricted low-frequency approximation. With the inclusion of gravitational radiation effects, the 
frequencies become complex- valued, thus potentially removing the singular structure of Kojima's equation. But even in the 
case of real- valued frequencies, it has been recently shown by Lockitch, Andersson & Friedman (2001) that for a non-barotropic 
uniform density model, in addition to the continuous spectrum there also exists a single mode solution with the eigenfrequency 
lying outside the continuous band. It is this mode that represents the relativistic r-mode for non-barotropic stars. 

In this paper, we extend the search of r-mode solutions to stars with various polytropic and realistic equations of state 
(EOS). We shall show that in addition to the continuous part, the eigenvalue equation always admits a single mode solution, 
at least for I — 2. However, for some stellar models, depending on the polytropic index n and on the compactness, the 
frequency of this solution lies inside the continuous band and is associated with a divergence in the fluid perturbation at 
the singular point. This is clearly not acceptable, and therefore we have to discard such solutions as being unphysical. As 
a logical consequence, we conclude that in those cases, there do not exist any r-modes, at least within the low-frequency 
approximation. In an independent work, Yoshida (2001) has come to a similar conclusion. He showed that even when studied 
in the post-Newtonian approximation, some polytropic models do not admit any r-modes. 

For realistic EOS, the existence of r-modes depends on the average polytropic index of the high-density regime. For very 
stiff EOS, the neutron star models can exhibit r-modes throughout the complete physically acceptable mass range, whereas 
for the very soft EOS, none of the neutron star models does. For EOS in the intermediate range, the existence of r-modes 
depends on the compactness of the stellar model. In addition to the mode calculations, we also use the time depend form of 
the equations. For those cases where we can find a physical r-mode solution, the Fourier spectrum of the time evolution does 
indeed show a peak at the appropriate frequency, whereas for those cases where we only have the unphysical mode, it does 
not. 

2 MATHEMATICAL FORMULATION 

Assuming that the star is slowly rotating with a uniform angular velocity fi, we neglect all terms of order higher than £1. In 
this approximation, the star remains spherical, because the deformation due to centrifugal forces is of order Q, 2 . Thus, the 
metric can be written in the form 

* Following Lockitch et al. (2001) we call a stellar model barotropic if the unperturbed configuration obeys the same one-parameter 
equation of state as the perturbed configuration. 




(1) 




(2) 
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dsl = -e 2v dt 2 + e 2X dr 2 + r 2 (d6 2 + sin 2 9dcj> 2 ) - 2cur 2 sin 2 9dtd<j) , (3) 

where v, A and the "frame dragging" u> are functions of the radial coordinate r only. With the neutron star matter described 
by a perfect fluid with pressure p, energy density e, and four-velocity 

(u\ u r , u\ it*) = (e - ", 0, 0, fie - ") , (4) 

the Einstein equations, together with a one-parameter equation of state p = p(e), yield the well-known TOV equations plus 
an extra equation for the function zu, defined as 

zu := Q, — u . (5) 
To linear order, this equation is 

zu" — ^4irre 2X (p + e) — — ^ zu' — 16Tve 2X [p + e) w = . (6) 
In the exterior, it reduces to 

zu" H — zu' = , (7) 
r 

for which we have the solution (Hartle 1967) 
2 T 

w = (8) 

with J being the total angular momentum of the neutron star. Equation (H) has to be integrated from the centre of the star 
to its surface R, where it has to match smoothly the exterior solution (H). With the angular momentum given by (Hartle 
1967, Glendenning 1997) 



8tt 



R 



J—~ir I "(p + tjzudr, (9) 

3 Jo 

the matching condition becomes 

R 4 zu'(R) = 6J , (10) 
and with Eq. (§) 

zu'(R) = ^-(tt-zu(R)) . (11) 
R 

If we focus on pure axial perturbations, the perturbed metric can be written in the following form: 
ds 2 = ds 2 ) + 2j2(hl™{t,r)dt + h[ m (t,r)dr) (-^^d6 + smOdeYimd^ , (12) 

l ,m 

where Yj m = Yi m (8, <j>) denote the scalar spherical harmonics. In addition, the axial component of the fluid velocity perturbation 
can be expanded as 

^{p + e)(8u e ,8u") =e"^l7 ,m (t ) r)(-^^ I 8infla e y !m ) . (13) 

l,m 

Einstein's field equations then reduce to four equations for the three variables /i m , h 1 ™ and U lm (Kojima 1992). 

As an alternative, we can use the ADM-formalism (Arnowitt, Deser & Misner 1962) to derive the evolution equations 
for the axial perturbations. The usefulness of this formalism for the perturbation equations of non-rotating neutron stars has 
been showed in Ruoff (2001). This formalism can be taken over to rotating stars and we can deduce equations describing the 
evolution of purely axial oscillations of slowly rotating neutron stars in terms of metric and extrinsic curvature variables. We 
should mention that our derivation starts with the complete set of perturbations, including both polar and axial perturbations. 
Only at the end do we neglect the coupling between the two parities and focus only on the axial equations. In the Regge- 
Wheeler gauge, there are just 2 non-vanishing axial metric perturbations and 2 axial extrinsic curvature components. In the 
notation of Ruoff (2001), they metric components are given by 



(fle,M = e"- A ^#r (-^^,sin0d fl F !m ) , (14) 
(hr$,hr4,) = e A -^l/ 4 im (-M|i,sin0d 9 r ;m ) , (15) 
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and the extrinsic curvature components read 

(k T e,k^) = \e x Y, K ^(~^f^ dd ^) > (18) 



2 " V sin6> 

l.m 

kee k e <t, \ 1 -a r ^(m . fl ( -sin6>~ 2 X !m W im \ 

k<f>e k^ I 2 ^ \ Wim Mm I 

l,m 

where Wi m and Xi m are abbreviations for 

Wlm = (flg + i(* + 1)) ^™ , (18) 

X ;m = 2{de-cot9)d 4> Y lm . (19) 
For the fluid velocity perturbation, we use the covariant form 

{5u e ,5u^) = e~ A y^ u l a m (t,r) ( t-^t , sinOdeYim) . (20) 

l,m 

The relation to the expansions ( |l^ ) and ( |l3| is given by (from now on we omit the indices I and m): 

h = e"- x K 6 , (21) 
fci = e A ~^4, (22) 
U = 4ire- x - v ( P + e)(u :i - K 6 ) . (23) 

We obtain the following quite simple set of evolution equations for the variables T4, K3, K§ and 113: 



2y-2A 



9 \ „ imr 2 



(24) 



( d ^- \ 1(1 + 1) -2 ^ 2im , _ 2A 



, w , imc;J J< 6 = Vl - jy—r [lu'K 3 - 167r(fi - w)(p + e)«s] , (26) 

(^+-fi)«a = ' («.-*,). (27) 

Furthermore, we have one momentum constraint: 

lfcr(p + £ ) W3 = K + ltf 3 _ ^ + 1 2 )- 2 Jf 6 - |^e--y 4 . (28) 
r t(t + 1) 

These equations are completely equivalent to the axial parts of Eqs. (20), (24), (25) and (27) of Kojima (1992) when the 
coupling to the polar equations is neglected. 



3 LOW-FREQUENCY APPROXIMATION 

The above evolution equations should not only describe the r-modes but also the axial ui-modes, which have much higher 
oscillation frequencies. If we want to focus on the r-modes only, we can use the fact that from Eq. ([!]), it follows that the 
r-mode frequency a is proportional to the star's angular velocity fi. Hence, we can require that in our evolution equations 
([m]) - ([27]), the time derivative dt be proportional to the r-mode frequency a or, equivalently, to £1 (Kojima 1997, 1998). In 
this case, we can order the perturbation variables in powers of Q (Lockitch et al. 2001) as 

u 3l K 3 ,K 6 ~ O(l), 

Va ~ O(fi) • (29) 

Keeping only the lowest order terms, we can neglect terms proportional to V4 in the evolution equation (^) and in the 
constraint which then read 

K' c> + (V - A' - ^ Ke - e 2X K 3 = , (30) 

K's + ^Ks- l{l+ ^~ 2 K 6 = 16n(p + e)u 3 . (31) 

These can be easily combined to give a single second-order differential equation for Kg. However, it is more convenient to 
write this equation in terms of the variable Hq — e v ~ x K%: 
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e h'a — Airr(p + e)h' + 



8n(p + e) 



AM 

r-3 



1(1+1) 



16ne v A (p + e)it3 



together with the evolution equation for M3 



d 



+ 



2vj 



1(1 + 1) 



- v h 



1*3 



(32) 



(33) 



At this point, it is worth making some comments on this approximation. The full set of axial equations (|24|) - (^) is a 
hyperbolic system describing the propagation of gravitational waves, which are excited on the one hand by the fluid motion 
(r-modes) and on the other hand by the curvature of spacetime itself (ui-modes). With the above approximation, we have 
completely suppressed the wave propagation, and the resulting equations now correspond to a Newtonian-like picture, where 
the fluid oscillations are acting as a source in the equation of the gravitational field. As now being described by a Poisson-like 
elliptic equation, implying an infinite propagation speed, the gravitational field ho reacts instantaneously on any changes 
in the source M3. Of course, this picture is only an analogy, since the metric variable ho corresponds to a post-Newtonian 
correction of the gravitational field and vanishes completely in the Newtonian limit. 

Furthermore our derivation of this approximation is only valid for non-barotropic stars. This is because in general we 
cannot start from decoupling the polar and axial equations in the first step as we did. Instead, when we apply the low- 
frequency approximation, we actually have to start from the full coupled system of equations, including both polar and axial 
perturbations. If we then do the same ordering in powers of fi, we also would have some polar variables of order 0(1), namely 
the remaining two fluid velocity components, coming from 5u r and the polar part of the angular components (Sue, 8u$), and 
the (rt) component of the metric, usually denoted by Hi. It turns out that the polar constraint equations can be combined 
to give a single constraint for Hi, which can be reduced to 



(r-ri)Hi=o, 



with 



p + € dp 



(34) 



(35) 



p de 

the adiabatic index corresponding to the unperturbed configuration and Fi the adiabatic index of the perturbed configuration, 
which in general differs from F. This is the case for non-barotropic stellar models, and Eq. ( |34| ) can only be satisfied if Hi 
vanishes. But this automatically implies that the polar fluid perturbations vanish, too, leaving thus only the axial equations, 
given above. In the barotropic case, it is F = Ti and the constraint for Hi is trivially satisfied, even for nonzero Hi . But this 
has as a consequence that the coupling between the polar and axial mode cannot be neglected, which means that there cannot 
exist pure axial mode solutions, since any kind of pure axial initial data will through the coupling automatically induce polar 
fluid oscillations. Hence, our analysis is strictly valid only for non-barotropic stellar models. 

As a further approximation, we could completely neglect all the metric perturbations. With this so-called relativistic 
Cowling approximation, we would be left with a single evolution equation for the fluid variable U3: 



d 

dl U3 



2zu 



1(1 + 1) 



u 3 



(36) 



From this equation we can immediately deduce that the various fluid layers are decoupled from each other, which means that 
each layer has its own real oscillation frequency given by 



a — —m 



Q- 



2vj 



1(1 + 1) 



(37) 



In the Newtonian limit (uj — » Q), this reduces to the well known result for the frequency of the r-mode given in Eq. (|l|). It 
should be pointed out that in the relativistic case, the presence of the frame dragging uj destroys the occurrence of a single 
mode frequency and gives rise to a continuous spectrum, at least to this order of approximation. Of course, it has been argued 
that this might be a pure artefact of the approximation, and the continuous spectrum may disappear as soon as certain 
approximations are relaxed. 

Let us therefore return to the low-frequency approximation, which is less restricted than the Cowling approximation, and 
see whether or not we can find real mode solutions. To this end, we assume our variables to have a harmonic time dependence 



u 3 (t,r) = u 3 (r)e lCTt , 
h (t, r) = ho(r)e~ lat 



(38) 
(39) 



Note that for the sake of notational simplicity, we do not explicitely distinguish between the time dependent and time 
independent form of the variables. With this ansatz, we assume the r-mode frequency a to be positive for positive values of 
m in contrast to the definitions in Eqs. (^) and ([!]). From Eq. ( ^j| ) we find that 
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2mvj 



«3 



2mw + 1(1 + l)(cr - mfi) 
which can be used to eliminate U3 in Eq. 
2mvj 



a — mQ, + 



1(1 + 1) 



— 2A ,11 

e h 



ho , 

_), yielding 
4:Tvr(p+ e)h' — ( 87r(p- 



4M + i; 



+ 167r(p + e)(a - mQ.)h Q = 



(40) 



(41) 



With appropriate boundary conditions, this equation represents an eigenvalue problem which should yield one, or possibly 
several, distinct eigenmodes. However, as was at first pointed out by Kojima (1997,1998), it might occur that the denominator 
in Eq. ( ^c| ) can become zero at some point inside the star, and the resulting eigenvalue problem becomes singular at this point. 
If the zero of the denominator lies outside the star, the eigenvalue problem is regular, since outside 113 = and Eq. ( |3^ ) can be 
directly solved without using Eq. ( pto| ) . The zeroes of the denominator occur if the frequency a lies in an interval determined 
by the values of vj at the centre and the surface, which we denote by vj c and vj s , respectively: 



ml! 1- 



2vj s 



£11(1 + 1) 



< a < mil 1 



2vj c 



ni(i + i) 



(42) 



By comparison with similar results from fluid dynamics, Kojima concluded that this equation should have a continuous 
spectrum with the frequency range given by ([l2|). This was put on a rigorous mathematical footing by Beyer & Kokkotas 
(1999). However, they could not exclude that there might not exist additional isolated eigenvalues, which would correspond 
to true mode solutions. 

We will now show that there actually exist such solutions, even though in some cases they are unphysical since the 
corresponding fluid perturbations would be divergent at the singular point. To make things look simpler, we can rescale 
Eq. ( J4l| ) and make it independent of Q and m. Following Lockitch et al. (2001), we introduce a normalized frequency 

(43) 



a = -1(1 + 1) (i--^tt 
2 y ' \ raS] 



and rewrite Eq. ([y]) as 

(q — vj 



1 , . s 4M 1(1 + 1) 
4nr(p + e)h' - ( 8tt(p + e) - — + v 2 ; 



/to 



+ 16tt(p + e)ah = 



where 

vj := zu/n . 

Equation (UoJ) then reads 



M3 = 



-e ho 



(44) 



(45) 



(46) 



vj — a 

Equation (^) becomes singular if a lies in the interval limited by the values of vj at the centre and at the surface of the star, 
i.e. if 



vj c < a < vj„ 



(47) 



For a solution to be acceptable, it must be regular at the origin, which amounts to ho(0) = 0, and it must vanish at infinity. 
As already mentioned above the integration of Eq. ([ll]) is straightforward, if the singular point lies outside the star. It is only 
when the singular point lies inside the star that some care has to be taken. 

Let us now assume the singular point r = ro lie inside the star. An analysis of Eq. ( [ll| ) (Andersson 1998b) shows that 
the singular point is a regular singularity, which admits a Frobenius expansion of the form 



+ . 



h (r-r ) = A (ai(r - r ) + a 2 (r - r ) 

B [(ai(r - r ) + a 2 (r - r ) 2 + . . .) In \r - r \ + b + 61 (r - r ) + b 2 (r - r ) 



rof + 



(48) 



Even though the solution is finite and smooth at the singular point r — ro, its derivative diverges, because of the logarithmic 
term. Moreover, if we want to compute the associated velocity component U3, we find that unless ho(r — ro) = 0, it will blow 
up at the singular point r = ro- Therefore, we conclude that the coefficient B has to vanish in order to obtain a physical 
solution and we are only left with the first power series starting with the linear term a±(r — ro). This yields vanishing ho at 
r = ro and therefore 113 can be finite at this point. The question is whether there are solutions satisfying both ho(r — ro) = 
and the appropriate boundary conditions at the centre and at infinity. We will now show that this cannot be the case. 
Suppose that ho(r = ro) = and h' (r = ro) > 0. We can then integrate Eq. ( fill ) from ro to r > ro: 



h'o(r) = h' (r ) + J <~ 

'r 



4nr(p + e)h' + ( 8ty( p + e) 3- 



4M 1(1 + 1) 



h + 16tt(p + e)- 



-hc 



dr. 



(49) 
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Since it is zu — a > for r > ro, each coefficient in the integral is strictly positive, hence we will have h' (r) > for all r > ro; 
i.e. ho keeps increasing asr-t oo, which is clearly incompatible with our requirement that ho vanish at infinity. Of course, the 
same argument holds if h' < at r = ro, with ho keeping decreasing. Hence, it follows that ho{r = ro) 7^ 0, but this means 
that we cannot have a vanishing coefficient B. Therefore, our solution will always contain the logarithmic term, which means 
that the associated velocity perturbation 113 is divergent at this point. This is clearly unphysical. We thus conclude that it 
is in principle possible to find eigensolutions to Eq. ([lij), however, if the associated eigenfrequencies lie inside the continuous 
band, the solutions become singular and have to be discarded on physical grounds. 

ft has been shown by Lockitch et al. (200f ) that for the existence of mode solutions, the allowed range of the eigenvalues 
a is bounded from below by vj c : 

vj c < a < 1 . (50) 

However, based on our above argumentation, we can further restrict this interval for the physically allowed eigenmodes to 
have as lower limit the value of w a '- 

uj s < a < 1 . (51) 



4 NUMERICAL RESULTS 

The numerical integration of Eq. (^) can be easily accomplished if the singular point lies outside the star, since in the exterior 
U3 = and we can therefore use the non-singular Eq. ( |32[ ) for the integration toward the outer boundary. If the singular point 
lies inside the star, we initiate our integration with a regular solution at the origin and integrate outward close to the singular 
point, where we match the solution to the expansion (^8|), i.e. we compute the leading coefficients 60 and ai. This gives us 
the new starting values to the right of the singular point and we can continue the integration up to a finite point outside the 
star, where we test if the solution satisfies the correct boundary condition. We mention again that for the integration outside 
the star we take the non-singular Eq. (p^), with U3 set to zero. 

We have performed mode calculations for sequences of uniform density and polytropic stars. For the uniform density 
models, the eigenfrequency a always lies outside the range of the continuous spectrum and therefore the associated eigen- 
functions do not exhibit any singularities. In Fig. we show the normalized eigenvalues a for I = 2 as a function of the 
compactness M/R together with vj c and w s , marking the boundaries of the continuous spectrum. For larger I (not shown), 
the eigenvalues a decrease and converge to cfr s , but stay always above -6j s . ^ In Fig. [], we show the eigenfunctions ho and 113 
for a uniform density model with a compactness of M/R = 0.153 and corresponding mode frequency a — 0.89806. Close to 
the centre of the star, the fluid perturbation 1*3 is proportional to r l+1 , but as it approaches the stellar surface, it grows much 
stronger, which comes from the denominator in Eq. (W6J) becoming very small. 



4.1 Polytropic models 

For polytropic models, obeying an equation of state of the form 

p = Ke 1+1/n (52) 

with polytropic index n, we obtain a quite different picture. For a polytropic index n = 1, as it is for instance shown in Fig. ^, 
it is only for the less compact stellar models that the eigenfrequencies a lie outside the continuous spectrum and therefore 
represent physical mode solutions. However, they are already that close to the upper boundary of the continuous spectrum 
■uj s that in Fig. ^ they cannot be distinguished any more. For more compact models, the eigenfrequency eventually moves 
inside the domain of the continuous spectrum, which means that the singular point now lies inside the star. This happens 
for a compactness of about M/R « 0.085. As discussed above, at the singular point the mode solution for ho exhibits an 
infinite slope and the associated fluid perturbation U3 diverges. Therefore, we have to discard them as being unphysical mode 
solutions. In Table [jj we have listed some polytropic models with their physical parameters and the eigenvalues a for I = 2 and 
/ = 3. The frequencies which are marked by an asterisk lie inside the continuous band and therefore correspond to unphysical 
mode solutions. For 1 = 2, only models I and 2 permit physical modes, whereas for I = 3, the modes are unphysical for all 
the stellar models. We should also mention that all our values are in perfect agreement with those previously obtained by 
Andersson (1998b). 



t Note that in Table 1 of Lockitch et al. (2001), there is a systematic error in their given values of a, which are too large by about 
5 per cent. This might be a consequence of a misprint in their Eqs. (5.2), (5.4) and (5.7), where the terms (1 — 2Mo/R) 1 ^ 2 and 
(l - 2M /R(r/R) 2 y /2 got confused. 
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Table 1. List of polytropic stellar models with n = 1 and k = 100 km 2 . 



Model 


e [g/cm 3 ] 


M [M Q ] 


R [km] 


M/R 






a(l = 2) 


a(l = 3) 


1 


1.0 x 10 14 


0.120 


12.32 


0.014 


0.96168 


0.99237 


0.99254 


0.98523* 


2 


5.0 x 10 14 


0.495 


11.58 


0.063 


0.83048 


0.96431 


0.96453 


0.92446* 


3 


1.0 x 10 15 


0.802 


10.81 


0.109 


0.70420 


0.93407 


0.93362* 


0.84895* 


1 


5.0 x 10 15 


1.348 


7.787 


0.256 


0.28377 


0.80236 


0.72579* 


0.44782* 


5 


1.0 x 10 16 


1.300 


6.466 


0.297 


0.14214 


0.74960 


0.52932* 


0.25301* 



To assess how the existence of a physical mode solution depends on the polytropic index n, we have computed modes 
for stellar models with fixed compactness M/R but with different values of n, ranging from to 1.5. The results are depicted 
in Fig. ^, where we show a as a function of n for / = 2. For small values of n, i.e. for stiff equations of state, the mode 
eigenfrequency lies outside the range of the continuous spectrum. But as n is increased, what corresponds to softening the 
equation of state, the mode frequency eventually crosses the boundary and migrates inside the continuous spectrum. This 
happens at n ~ 0.8, but for larger values of the transition point moves to smaller values of n. Actually, it is not the mode 
frequency a which moves towards the boundary of the continuous spectrum zu a as n is increased, it is rather the boundaries 
of the continuous spectrum which start to expand, and zb s approaches the mode frequency a, which more or less hovers at a 
constant value. For n — 0, the uniform density models, the range of the continuous spectrum (the shaded area in Fig. ^) is 
the smallest, and probably it is only in this case that one can find eigenvalues for quite large I, if not for all We should also 
mention that for each polytropic model, there seems to exist a maximal value of I, beyond which there are neither physical 
nor unphysical mode solutions. The frequency a of the unphysical mode solution quickly approaches vj c as I is continuously 
increased. For I larger than the critical value, where a(l) = w a , we could not find any mode solution at all. For the n = 1 
models of Table 1, this happens already for I = 4. 

To check and corroborate our above results, we also numerically evolved the time dependent equations ( p2] ) and (^) and 
took Fourier transforms of the resulting evolution. For initial data representing the physical mode solution of Fig. ^, the time 
evolution indeed gives a single frequency signal at each point inside and outside the star. In this case, there is no sign of a 
continuous spectrum at all, and all the fluid layers oscillate in a uniform manner. This is shown in Fig. [s| where for both the 
fluid variable «3 (left panel) and the metric variable ho (right panel), there is one single peak, which is independent of the 
location r. 

If we now choose arbitrary initial data, as for instance in Fig. ^, we expect the power spectrum at a given location r to 
consist of two peaks: One which is independent of the location inside the star and represents the eigenmode, and another peak 
which varies between the boundaries determined by ro c and vj s as one moves throughout the star. This is how the continuous 
spectrum should show up in the Fourier transform. For the fluid variable «3, the power spectrum of the evolution indeed 
confirms our expectations, as is shown in the left panel of Fig. [?]. 

However, the spectra of ho show that for locations closer to the stellar surface, the peaks corresponding to the continuous 
spectrum are smaller by several orders of magnitudes compared to the peak representing the eigenmode. For 113, the peaks 
are of the same order of magnitude. And outside the star, ho shows only the eigenmode peak, and no sign of the presence of 
the continuous spectrum, which should reveal itself as a superposition of all the frequencies in the range between ro c and zu s - 
In therefore seems to be invisible for an external observer. We should note that those spectra are taken after a certain initial 
time, in which the system adjusts itself. If we had taken the Fourier transform right from t = 0, we would have obtained a 
clear sign of the continuous spectrum. 

Let us now turn or attention to the polytropic cases, where we can have unphysical mode solutions. We will present 
evolution runs for the stellar models 1 and 5 from Table |l| with 1 — 2. For model 1, the singular point lies outside the stellar 
surface and therefore there exists a physical mode solution. For model 5, the singular point lies inside the star, hence a < zu a 
and the associated eigensolution is unphysical. It should be noted that this model is also unstable with respect to radial 
collapse. 

For model 1, the physical mode solution can be used as initial data. As for the uniform density case, the numerical evolution 
of such data yields a purely sinusoidal oscillation with the expected r-mode frequency a. Therefore, the corresponding power 
spectrum is similar to Fig. For arbitrary initial data, we obtain a picture similar to Fig. ^. Note, that the values of a and 
vj s differ only by about 0.01 per cent. Still, with a high resolution run we can numerically distinguish these values, as is shown 
in Fig. |i| Here, we plot the power spectra of ho(t) and 113 (t) taken at the stellar surface. The spectrum of ho shows one single 
peak at the eigenfrequency a, whereas M3 shows two peaks at a and zb s - 

For model 5, things are quite different. Here, we cannot evolve initial data representing the unphysical mode solutions 
because the fluid perturbation would diverge at the singular point. Yet, if this solution still had some physical relevance, then 
arbitrary initial data should be able to excite this mode, and the power spectrum of the time evolution should show a peak at 
the corresponding frequency. However, this is clearly not the case, as can be seen in Fig. ^, where we show the late time power 
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spectra of the time evolution of 1*3 (left panel) and ho (right panel) for model 5. For the fluid variable 113, there is always one 
single peak, which varies for different locations r between the boundaries vjc and vj s . There is not even the slightest trace of 
a common peak at the expected value of a — 0.52932. 

For the metric variable ho (right panel), we essentially observe the same picture. Here, too, no common peak can be found 
at the expected mode frequency a, but curiously there is nevertheless an additional common peak for all locations with its 
frequency given exactly by m 3 . However, this peak does not show up in the power spectrum of 113, except, of course, directly 
at the surface. 

It is obvious that it cannot be a mode solution, since first of all the eigenvalue code does not give a solution for this 
particular frequency zb s or even in the close vicinity. Moreover, the time evolution shows a quite different behaviour compared 
to the case where a physical r-mode exists. In Fig. [n], we plot the time evolutions of ho outside the star for models 1 and 5 
of Table |l|. For model 1, where we have a physical r-mode, after some initial time the amplitude remains constant, whereas 
for model 5 the amplitude keeps decreasing with time and in this case fits very well a power law with an exponent of -2. For 
model 1, the dominant oscillation frequency is the corresponding r-mode frequency a, whereas for model 5 it is given by zu 3 . 
In both cases, the amplitude of the fluid perturbation 113 remains constant after some initial time. It now becomes clear, why 
we cannot observe the common peak at ro 3 in the fluid spectrum. The spectra are taken at late times, where the amplitude 
of ho and therewith its influence on U3 has considerably decreased. If we had taken the spectra at earlier times, we could 
observe a similar peak in the fluid spectrum, as well. 

We have no clear explanation what causes this additional peak, but we suppose that it comes from the behaviour of 
the energy density e at the surface. The peak is much more pronounced for polytropes with n < 1, since there the energy 
density e has an infinite slope at the surface. For n = 1, the slope e is finite and for n > 1, it is zero. In the latter case, the 
peak is strongly suppressed. Even for uniform density models, one can observe this additional peak, arising because of the 
discontinuity of the energy density at the surface. However, this peak is several orders of magnitudes smaller than the peak 
corresponding to the eigenmode, which is always present for uniform density models, and therefore hard to detect. 

It is instructive to compare the evolution of the same initial fluid perturbation for a uniform density model and a 
polytropic model, having the same compactness but without the latter admitting a mode solution. Since in the low-frequency 
approximation, there is no radiation which can dissipate the energy of the fluid, the total energy of the system should be 
conserved. However, we have observed that in the polytropic case the amplitude of ho is constantly decreasing, hence its initial 
energy has to be transferred to the fluid, whereas in the uniform density model, the energy should be shared between 113 and 
ho- This is indeed, what can be observed. In the uniform density case, the fluid amplitude does not change too much, but in 
the polytropic case, it shows a quite strong initial growth, accompanied by the strong decrease of ho- 

4.2 Realistic Equations of State 

Having found that for a quite large range of polytropic stellar models, there do not exist any physical r-modes, an obvious 
question is, whether or not realistic equations of state do admit physical mode solutions. To give an answer, we have investigated 
the collection of realistic equations of state which have been studied by Kokkotas & Ruoff (2001) for the radial modes. The 
relevant notations, references and data of the stellar models can be found in there. 

The results are quite unexpected and seemingly contradictory. When trying to compute the modes through Eq. (jiil), we 
find that for all the equations of state, the frequencies always lie inside the continuous band. Based on our above discussion, 
we therefore would have to discard them as being unphysical. It then would seem that none of the existing realistic EOS 
admits an r-mode, at least in the physically relevant range from about one solar mass up to the stability limit of each EOS. 
The surprise is now that the time evolution does show a quite different picture. Only for the EOS B (Pandharipande, 1971), 
G (Canuto & Chitre, 1974) and MPA (Wu et al., 1991) do the evolutions meet our expectations and show the continuous 
decrease in the amplitude of the metric variable ho, in accordance with the polytropic cases without r-modes. However, for 
all other EOS, the amplitude remains constant after a while, indicating that there is indeed a mode present. Only when 
approaching their respective stability limit do some EOS show a decay of the amplitude of ho- When obtaining the frequency 
through Fourier transformation, we find that it always lies inside the continuous spectrum, however, it does not coincide with 
the frequency found from the mode calculation. Instead, the frequency is in all cases given by the value of w close to the 
neutron drip point. 

How can we explain this discrepancy with our previous considerations? First, we would like to stress that it is not a 
numerical artefact of the time evolution, for convergence tests corroborate the presence of this mode. When examining the 
different EOS, we find that the EOS B, G and MPA are the softest ones, with a maximal polytropic index in the high density 
regime around n = 0.8. All others have indices less than 0.8, going down to n ~ 0.5 for the EOS I (Cohen et al. 1970) and 
L (Pandharipande et al. 1976). From Fig. ^, it becomes clear that it is just for polytropic models with n > 0.8 that the 
eigenvalue migrates inside the continuous band, and the r-mode therefore ceases to exist. For models with n < 0.8, we usually 
can find a physical r-mode, but this depends on the compactness of the model under consideration. 

At a density of 10 14 g/cm 3 , the effective polytropic index of any EOS is given by n ~ 2. As the EOS is approaching the 
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neutron drip point at a density of e f» 4 x 10 11 g/cm 3 , the EOS becomes softer and softer, i.e. n increases even further. Only 
for densities below the neutron drip point does the EOS stiffen again. This structure of the EOS is responsible of putting a 
low density layer (the crust) around the high density core. However, because of its low density compared to the core, this 
layer practically does not contribute to the total mass, its only effect is to slightly increase the radius of the neutron star. 
Thus, whether or not the EOS admits a mode, should be determined solely by the core. To assess this proposition, we can 
do the following. We replace the whole part of the EOS below 10 14 g/cm 3 by a smooth polytropic EOS with polytropic index 
of n = 2. By doing so, we obtain a model with practically the same mass, but a somewhat smaller radius. Computing the 
modes of the thus modified model, we find that the r-mode frequencies a actually do lie outside the continuous spectrum, if 
the average polytropic index of the core is less than 0.8, and the model is not too compact. However, a lies extremely close 
to the value of vj s . For softer EOS with an average index of n « 0.8, we still would not be able to find any physical modes. If 
we now go back and restore the outer layer, the r-mode frequency should not significantly change, because of the negligible 
gravitational influence of this outer layer. The only effect is the slight increase of the stellar radius R. But with R becoming 
larger, the value of vo a also increases and actually becomes larger than the r-mode frequency a, which then lies inside the 
continuous band. This is shown in Fig. [ll], where we plot the two density profiles for a stellar model based on the EOS WFF 
(Wiringa, Ficks & Fabrocini, 1988) and the same model with the low-density regime replaced by a n = 2 polytropic fit. For 
the polytropic fit, the zero of vj — a lies right outside the star, whereas for the complete realistic model, it is inside the star. 
In the latter case, however, the mode still exists, but it cannot have a purely harmonic time dependence any more. If this 
were the case, it had to be a physical solution of Eq. (^) , but it is clearly not since the frequency lies inside the continuous 
spectrum. We therefore conclude that this stable oscillation that can be seen in the time evolution is always a mixture of a 
mode and the continuous spectrum. We should mention that in our treatment of realistic EOS, we have assumed that the 
neutron stars consists entirely of a perfect fluid, even in the outer layer. This is certainly not true, instead a neutron star 
should have a solid crust, which clearly will modify the above results. However, this is beyond the scope of this work. 



5 SUMMARY 

We have performed both mode calculations and time evolutions of the pure axial perturbation equations for slowly rotating 
stars in the low-frequency approximation. Although the time independent equation (^) represents a singular eigenvalue prob- 
lem, admitting a continuous spectrum, it is nevertheless possible to find discrete mode solutions, representing the relativistic 
r-modes. If the mode frequency lies outside the continuous spectrum, the eigenvalue problem becomes regular, and the associ- 
ated solution represents a physically valid r-mode solution. If the eigenvalue lies inside the continuous band, the eigenfunction 
exhibits an infinite slope at the singular point, which is due to the presence of a r log |r | term in the series expansion. Moreover, 
the corresponding fluid perturbation U3 diverges at the singular point. Therefore, we conclude that these mode solutions are 
unphysical, and the only physically valid mode solutions are the ones where the associated frequencies a lie outside the range 
of the continuous spectrum. 

We have performed mode calculations for uniform density models, for various polytropic models and also for a set of 
realistic equations of state. In agreement with the results of Lockitch et al. (2001), we find that uniform density models 
generally admit r-modes for any compactness. For polytropic equations of state, however, the existence of physical r-mode 
solutions depends strongly on both the polytropic index n and the compactness M/R of the stellar model. The general picture 
is that the smaller n is, which corresponds to a stiffer equation of state, the larger is the compactness range where one can 
find physical mode solutions. For a given polytropic index n, one usually finds physical mode solutions for models with a small 
M/R ratio. As the compactness is increased, i.e. as the models become more relativistic, the mode frequency a decreases 
and starts approaching ro s . Eventually it crosses this point and migrates inside the range of the continuous spectrum, thus 
becoming unphysical, and no r-mode exists any more. 

When considered as a function of I, the r-mode frequency a is monotonically decreasing. For a uniform density model, a 
approaches vj s as I is increased. In polytropic models, this has the effect that it is even harder to find physical mode solutions 
for higher values of I, since a will much sooner cross the border vj a of the continuous spectrum. If I is large enough, it seems 
that the eigenvalue equation (Q) does not admit any mode solution at all, not even a singular one. 

We have verified our results by explicitly integrating the time dependent equations. The time evolutions for the models 
admitting an r-mode can clearly be distinguished from those without a discrete mode. In the former case, both the fluid 
perturbation 113 and the metric perturbation ho oscillate with a constant amplitude after some initial time. In the latter 
case, the amplitude of ho constantly decreases. The fluid amplitude, however, still remains at a constant level. This can be 
explained by a decoherence effect in the fluid oscillations, since the frame dragging causes each fluid layer to oscillate with a 
different frequency. Thus the initially uniform fluid profile becomes more and more disturbed because the fluid layers get out 
of phase, resulting in a continuously weakening of the strength of the fluid source term in Eq. ([52]) . When a physical r-mode 
exists, the system can oscillate in a coherent manner. 

When turning to realistic equations of state, the mode calculations yielded only frequencies lying inside the continuous 
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band, therefore being apparently unphysical. However, for most EOS the numerical time evolutions revealed the presence of a 
mode, but with the frequency still lying inside the continuous band and corresponding approximately to the value of vj at the 
neutron drip point. We explained this seemingly contradictory behaviour by making the core responsible for the existence of 
the r-mode. For, if we remove the outer layer, which does not have any significant gravitational contribution, we can indeed 
find eigenvalues a which lie outside the continuous band. However, they are very close to the upper limit of the continuous 
band w g . By adding the additional layer, we increase zb 3 such that it now becomes larger than a, which remains basically 
unaffected. Although now lying inside the continuous spectrum, the mode still exists, but it will be always associated with 
an excitation of the continuous spectrum. Most of the realistic EOS do admit r-modes in a certain mass range, but only the 
stiffest ones admit modes throughout the whole mass range up to the stability limit. The less stiff ones have a maximal mass 
model above which there are no r-modes any more, and for the softest EOS, namely EOS B, G and MPA, there are no r-mode 
for the whole physically relevant mass range. 

It should be kept in mind that all our results concerning the r-modes are obtained within the low-frequency approximation. 
It would be clearly much too early to infer any statements about the existence or non-existence of the r-modes in rapidly 
rotating neutron stars. And if the true EOS of neutron stars is rather stiff, and therefore would already admit r-modes within 
the low-frequency approximation, then the whole discussion about the singular structure would be irrelevant. But as the true 
EOS is not known yet, we cannot exclude it to be rather soft, and the appearance of the singular points has to be taken much 
more seriously. Still, it could still be seen as a mere artefact of the low-frequency approximation. But the work of Kojima & 
Hosonuma (1999) indicates that the inclusion of second-order terms in Q even increases the range of the continuous spectrum, 
which is responsible for the disappearance of the r-mode. They worked only in the Cowling approximation, but whether or 
not the inclusion of more higher order terms and the complete radiation reaction can restore the existence of the r-modes in 
all cases is still an open question and deserves further investigation. As a next step in this direction, we will investigate in a 
subsequent paper the full set of axial equations (Eqs. (H) - (I27])), which contains the radiation reaction. 
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M/R 



Figure 1. Boundaries of the continuous spectrum tb c and ro s together with the r-mode frequency a as a function of compactness M/R 
for uniform density models. The mode frequency a lies always outside the continuous band (shaded area). 




Figure 2. Eigenfunctions ho and us for a uniform density model with compactness M/R = 0.153. The corresponding r-mode frequency 
is a = 0.89806. 
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Figure 3. Boundaries of the continuous spectrum ro c and w s together with the r-mode frequency a as a function of compactness M/R 
for polytropic n = 1 models. For M/R < 0.085, the mode frequency a lies outside the continuous band (shaded area), but differs from 
the value of ro s not more than 0.01 per cent. For M/R > 0.085 it migrates inside the band. 




Figure 4. The boundaries of the continuous spectrum vj c and ■w s together with the r-mode frequency a as a function of the polytropic 
index n for stellar models with compactness of M/R = 0.2. The case n = corresponds to a uniform density model. 
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Figure 5. Power spectrum of the time evolution of a mode solution (initial data of Fig. |2j) for the uniform density model. For both the 
fluid us (left graph) and the metric variable ho (right graph), the spectrum shows a single peak at the expected frequency of a = 0.89806. 
In each graph, the power spectrum has been taken at 5 different locations inside the star. For ho, there is an additional one from outside 
the star. The spectra have been individually rescaled in order to give clearer graphs. 




Figure 6. Arbitrary initial data for the uniform density model. 
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Figure 7. Late time power spectrum of the time evolution of arbitrary initial data (initial data of Fig. for the uniform density model. 
At each location r inside the star, the fluid 113 (left graph) shows two peaks, one corresponding to the r-mode frequency a = 0.89806 
and another corresponding to the value of vj at this particular location r. The spectrum of the metric variable /io (right graph) is quite 
similar, however, for the spectra corresponding to locations close to the surface and outside the star, the influence of the continuous 
spectrum becomes negligible and the mode dominates. Again, the amplitudes of the spectra are arbitrarily rescaled. 




Figure 8. Power spectrum of the time evolution of arbitrary initial data for the polytropic model f , taken at the surface. Both U3 and 
ho show a peak at the expected r-mode frequency a = 0.99254. But only 113 shows an additional peak at m = vj s = 0.99237. 
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Figure 9. Late time power spectrum of the time evolution of arbitrary initial data for the polytropic model 5. For this model, we do not 
expect a physical mode to exist. The frequency of the unphysical mode is given by a = 0.52932. It is clear that at this location, neither 
U3 nor ho show a peak in the power spectra. At each location, the spectrum of 113 shows a single peak corresponding to the respective 
value of ■uj. However, in addition to those variable peaks, the spectra of h§ reveal a common peak with the frequency given by vd s . This 
peak can be traced back to the kink in the energy density profile at the surface. 



e„ = 10 g/cm 




t in ms 



Figure 10. Logarithmic plot of the time evolution of ho outside the star for the polytropic models 1 and 5 of Table 1. For model 1, 
which admits a physical mode solution, the oscillation amplitude remains constant, and ho oscillates with the r-mode frequency a. Model 
5 does not admit a physical mode solution and here, the amplitude of ho decreases with the decay fitting very well a power law with an 
exponent of -2. The oscillation frequency is given by the values of ro s . 
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Figure 11. It is shown the density profiles near the surface for neutron star model obtained from EOS WFF with a central density of 
2 X 10 15 g/cm 3 . In the one case we include the low density part, whereas in the other case this part is replaced by a polytropic fit. Also 
shown is the location where vd = a. For the polytropic fit, it lies outside the star, therefore representing a physical mode, but for the 
complete realistic model, it would lie inside. 
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